# ---------------------------------------------------------------- #
# The Journal of Politics
# When Do Citizens Resist The Use of AI Algorithms in Public Policy?
#
# Replication Code for Theory Validation (Figure A-4, Online Append)
#
# Author: Shir Raviv
# Date: Feb 16, 2025
# ---------------------------------------------------------------- #
# Clear Environment and Free Memory -------------------------------------------
rm(list = ls())
gc()

# Load Required Packages --------------------------------------------------------
required_packages <- c(
  "dplyr", "readr", "expss", "scales", "ggeffects", "lme4", "forcats",
  "vtable", "xtable", "stargazer", "texreg", "broom", "rlang",'stringr',
  "rio", "ggplot2" )

# Install (if needed) and load each package
for (pkg in required_packages) {
  if (!require(pkg, character.only = TRUE)) {
    install.packages(pkg, dependencies = TRUE)
    library(pkg, character.only = TRUE)
  }
}


# Set Working Directory and Define Export Paths -------------------------------
#setwd("insert here.. /JOP Replication")
plotpath <- "figures/"

# Load Data ---------------------------------------------------------------------
df <- readRDS("theory_validation_data.RDS")

# Generate Outcome Variables ----------------------------------------------------
df <- df %>%
  mutate(
    # Outcome indicators for each type of response
    sacnt_col_resp = case_when(
      res_oucome == "Sanction collectives" ~ 1,
      res_oucome != "Sanction collectives" & !is.na(res_oucome) ~ 0
    ),
    sacnt_ind_resp = case_when(
      res_oucome == "Sanction individuals" ~ 1,
      res_oucome != "Sanction individuals" & !is.na(res_oucome) ~ 0
    ),
    assis_col_resp = case_when(
      res_oucome == "Assist collectives" ~ 1,
      res_oucome != "Assist collectives" & !is.na(res_oucome) ~ 0
    ),
    assis_ind_resp = case_when(
      res_oucome == "Assist individuals" ~ 1,
      res_oucome != "Assist individuals" & !is.na(res_oucome) ~ 0
    ),
    # Composite outcomes
    sacntioning_resp = case_when(
      res_oucome %in% c("Sanction collectives", "Sanction individuals") ~ 1,
      res_oucome %in% c("Assist collectives", "Assist individuals") ~ 0
    ),
    collectives_resp = case_when(
      res_oucome %in% c("Sanction collectives", "Assist collectives") ~ 1,
      res_oucome %in% c("Sanction individuals", "Assist individuals") ~ 0
    )
  )

# Run Regressions ---------------------------------------------------------------
# Regression for "collectives_resp"
lm_collectives <- lm(collectives_resp ~ factor(decision_label) + factor(workerId), data = df)
lm_collectives <- tidy(lm_collectives, conf.int = TRUE) %>%
  filter(!grepl("TP_", term))

# Regression for "sacntioning_resp"
lm_sacntioning <- lm(sacntioning_resp ~ factor(decision_label) + factor(workerId), data = df)
lm_sacntioning <- tidy(lm_sacntioning, conf.int = TRUE) %>%
  filter(!grepl("TP_", term))

# Combine regression results
coef_df <- bind_rows(
  lm_collectives,
  lm_sacntioning,
  .id = "type_outcome"
)

coef_df$type_outcome <- factor(coef_df$type_outcome,
                               levels = c(1, 2),
                               labels = c("Collectives", "Sanctioning"))


coef_df <- coef_df %>%
  filter(term != "(Intercept)") %>%
  mutate(term = str_remove(term, "factor\\(decision_label\\)"))

# Add Theoretical Classification
coef_df <- coef_df %>%
  mutate(theory_classification = case_when(
    term %in% c("study", "teacherprom", "visa", "foodstamps", 
                "covid_tests", "mental_service", "caseworker_coaching") ~ "Assist individuals",
    term %in% c("streetlighting", "community_centers", "educ_funding", 
                "shelters", "firestations") ~ "Assist collectives",
    term %in% c("lockdowns", "illegalwork", "drugtests", "police_patrols", 
                "illegal_construction") ~ "Sanction collectives",
    term %in% c("monitoring", "childabuse_investigation", "sentence", 
                "restraining", "fireteacher") ~ "Sanction individuals"
  ))
coef_df$theory_classification <- factor(coef_df$theory_classification)

# Create additional classification variables
coef_df <- coef_df %>%
  mutate(
    theory_assisting = case_when(
      theory_classification %in% c("Assist collectives", "Assist individuals") ~ "Assisting",
      theory_classification %in% c("Sanction collectives", "Sanction individuals") ~ "Sanctioning"
    ),
    theory_collectives = case_when(
      theory_classification %in% c("Assist collectives", "Sanction collectives") ~ "Collectives",
      theory_classification %in% c("Assist individuals", "Sanction individuals") ~ "Individuals"
    )
  )
coef_df$theory_assisting <- factor(coef_df$theory_assisting)
coef_df$theory_collectives <- factor(coef_df$theory_collectives)

# labels 
coef_df <- coef_df %>%
  mutate(label = case_when(
    term == "educ_funding" ~ "Deciding which schools should receive extra funding",
    term == "streetlighting" ~ "Deciding where to place street lighting",
    term == "firestations" ~ "Deciding where to place fire stations",
    term == "community_centers" ~ "Deciding where to open community resource centers",
    term == "shelters" ~ "Deciding where to build shelters",
    term == "visa" ~ "Determining visa immigration for an applicant",
    term == "teacherprom" ~ "Deciding which teachers to promote",
    term == "study" ~ "Deciding which pupils should receive study assistance",
    term == "foodstamps" ~ "Deciding which individuals should receive food stamps",
    term == "illegalwork" ~ "Deciding where to increase enforcement of illegal work",
    term == "illegal_construction" ~ "Deciding where to increase enforcement of illegal constructions",
    term == "police_patrols" ~ "Deciding where police forces should patrol",
    term == "lockdowns" ~ "Deciding in what areas to impose a lockdown",
    term == "drugtests" ~ "Deciding at which schools to conduct drug/alcohol tests",
    term == "restraining" ~ "Determining a restraining order",
    term == "sentence" ~ "Determining a sentence for a defendant",
    term == "fireteacher" ~ "Deciding which teachers to fire",
    term == "monitoring" ~ "Deciding which residents the police should monitor",
    term == "childabuse_investigation" ~ "Deciding which families to investigate for child abuse",
    term == "covid_tests" ~ "Deciding which individuals should get a COVID vaccination",
    term == "mental_service" ~ "Deciding which residents should receive mental health assistance",
    term == "caseworker_coaching" ~ "Deciding which families to provide caseworker coaching"
  ))

coef_df$theory_classification <- factor(coef_df$theory_classification,
                                        levels = c("Assist collectives",
                                                   "Assist individuals",
                                                   "Sanction collectives",
                                                   "Sanction individuals"))

# Figure  ------------
plot <- ggplot(coef_df, aes(y = fct_rev(label), x = estimate,
                            colour = type_outcome, fill = type_outcome)) +
  facet_grid(theory_classification ~ type_outcome, scales = "free_y", space = "free_y") +
  geom_vline(xintercept = 0, linetype = "dashed", color = "gray80") +
  geom_linerange(aes(xmin = conf.low, xmax = conf.high), size = 0.3, alpha = 0.8) +
  geom_linerange(aes(xmin = estimate - std.error * 1.645,
                     xmax = estimate + std.error * 1.645), size = 0.95, alpha = 0.8) +
  geom_point(size = 2) +
  scale_fill_manual(values = c("#397685", "#719482")) +
  scale_colour_manual(values = c("#397685", "#719482")) +
  theme_minimal() +
  theme(
    panel.grid.minor = element_blank(),
    panel.grid.major.y = element_blank(),
    legend.position = "none",
    strip.text.y = element_text(angle = 90, size = 9),
    strip.text.x = element_text(size = 10),
    strip.placement = "outside",
    axis.text.x = element_text(size = 9),
    axis.text.y = element_text(size = 9, hjust = 1),
    axis.title.x = element_text(size = 10, margin = margin(t = 10)),
    panel.spacing.x = unit(2, "lines"),
    strip.background = element_blank(),
    axis.ticks.y = element_blank()
  ) +
  scale_x_continuous(expand = c(0, 0)) +
  labs(x = "Coefficient estimate", y = NULL, title = NULL)
plot

# Save the Figure 
ggsave(plot,
       filename = paste0(plotpath, "figA4.tif"),
       device = "tiff",
       height = 6, width = 10, dpi = 600)


ggsave(plot,
       filename = paste0(plotpath, "figA4.png"),
       device = "png",
       height = 6, width = 10, dpi = 600)


